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Abstract 

Max-algebra models of tandem single-server queueing systems with 
both finite and infinite buffers are developed. The dynamics of each 
system is described by a linear vector state equation similar to those 
in the conventional linear systems theory, and it is determined by a 
transition matrix inherent in the system. The departure epochs of a 
customer from the queues are considered as state variables, whereas its 
service times are assumed to be system parameters. We show how tran- 
sition matrices may be calculated from the service times, and present 
the matrices associated with particular models. We also give a repre- 
sentation of system performance measures including the system time 
and the waiting time of customers, associated with the models. As 
an application, both serial and parallel simulation procedures are pre- 
sented, and their performance is outlined. 

Key- Words: max-algebra, tandem queues, dynamic state equation, 
performance measure, parallel simulation algorithm. 



1 Introduction 

In the analysis of queueing systems, algebra models arise naturally from the 
recursive equations of the Lindley type, which present a formalism widely 
used for the representation of dynamics of a variety of queueing systems 
models. There are the recursive equations designed to describe the G/G/m 
queue, closed and open tandem queueing systems which may have both 
infinite and finite buffers, and queueing networks with deterministic routing 
(see e.g., [H [21 |3l [H [5l E] ) . These equations, which allow the dynamics of 
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a queueing system to be represented in a convenient and unified way well 
suited to analytical treatments, also provide the basis for the development 
of efficient procedures of queueing system simulation [U El |3l [7] . 

Since recursive equations often involve only the operations of arithmetic 
addition and maximization, they offer the prospect of the representation of 
queueing system models in terms of the max-algebra theory [HI [9l [TO]. The 
implementation of max-algebra allows one to rewrite the recursive equations 
as linear scalar and vector algebraic equations [11] , which are actually almost 
identical to those in the conventional linear system theory. The benefits of 
the max-algebra approach in the analysis of queueing systems are twofold: 
first, it gives us the chance to exploit results of the conventional linear 
algebra, which have been reformulated, and are now available in the max- 
algebra. The classical results already reformulated and proved in the max- 
algebra include, in particular, the solution of the eigenvalue problem, the 
Cayley-Hamilton theorem, and Cramer's rule [HI [9} \T0[ I12j. 

Other benefits have a direct relationship to computational aspects of 
simulation. In fact, the algebraic models of queues lead to matrix-vector 
max-algebra multiplications as the basis of simulation procedures [2]. New 
possibilities then arise in queueing system simulation to employ efficient 
computational methods and algorithms available in numerical algebra, in- 
cluding those designed for implementation on parallel and vector processors. 

In this paper we develop max-algebra models of open and closed tan- 
dem single-server queueing systems which may have both infinite and finite 
buffers, and we give related representations of system performance mea- 
sures. We start with preliminary algebraic definitions in Section 2 which 
also includes a technical lemma underlying the development of models in 
later sections. In Section 3, the dynamics of open and closed systems with 
infinite buffers is described by a vector state equation which is determined 
by a transition matrix inherent in the system. The departure epochs of a 
customer from the queues are considered as state variables, whereas its ser- 
vice times are assumed to be system parameters. We show how transition 
matrices may be calculated from the service times, and present the matrices 
associated with certain particular models. 

Section 4 extends the dynamic equation to cover open tandem systems 
with finite buffers, which operate under both manufacturing and communi- 
cation blocking rules. The representations of system performance measures 
including the system times and the waiting times of customers, associated 
with the models are given in Section 5. In Section 6, we present serial and 
parallel simulation algorithms based on the algebraic models, and outline 
their performance. Finally, Section 7 gives conclusions. 
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2 Preliminary Algebraic Definitions and Results 



In this section we briefly outline basic facts about matrix max-algebra, which 
underlie the algebraic models and methods of queueing system simulation, 
presented in the subsequent sections. Further details concerning the max- 
algebra and its applications can be found in survey papers [12]. A 
thorough theoretical analysis of this algebra and related algebraic systems 
is given in O E] . 

We start with the max-algebra of real numbers, which is the system 
(K, ©, 0) , where R = R U {e} with e = -oo, and 

X (B y = max(x, y), x ® y = x + y 

for any x, y € R . 

It is easy to see that these new operations, namely addition © and 
product (8) , possess the following properties: 

Associativity: x © (y © z) = (x © y) © z, 

X (g) (y (g) z) = (x (g) y) (g) z; 
Commutativity: x©y = y©x, x y = y x; 

Distributivity: x (g) (y © z) = (x (g) y) © (x (g) z); 

Idempotency of Addition: x © x = x. 

With e = 0, we further have 

Null and Identity Elements: x©e = effix = x, x(g)e = e(g)x = x; 
Absorption Rule: x(ge = e(g)x = e. 

Clearly, in the max-algebra these properties allow ordinary algebraic manip- 
ulation of expressions involving the max-algebra operations to be performed 
under the usual conventions regarding brackets and precedence of (g over 
ffi. 

Note finally that in the max-algebra, for each x G R, there exists its 
multiplicative inverse x~^ such that x (g) x~^ = x~^ (g) x = e. This is the 
usual arithmetic inverse which satisfies, in particular, the evident condition 

(x (g) y)~^ = x^'^ (g) y^^ 

for all X, y G R. 

2.1 Max-algebra of Matrices 

The scalar max-algebra is extended to the max-algebra of matrices in the 
regular way. Specifically, for any square {n x n) -matrices A = (aij) and 
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B = (bij) with entries in R, the elements of the matrices C = A ® B and 
D = AiS> B are calculated as 



and 



k=l 



where J2® denotes the iterated operation ©, i = l,...,n; j = l,...,n. 
Similarly, the multiplication of a matrix by a scalar, as well as the operations 
of both matrix-vector multiplication and vector addition may be routinely 
defined. 

As in the scalar max-algebra, there are null and unit elements in the 
matrix algebra, defined respectively as 





(e . 
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One can easily see that for any square matrix A, it holds 

S0A = A^S = £, E0A = A(S)E = A, S®A = A®£ = A. 

It is not difficult to verify that the other properties of scalar operations 
© and 0, with the exception of the commutativity of multiplication, are 
also extended to the matrix algebra. Similar to the conventional matrix 
algebra, matrix multiplication in the max-algebra is not commutative in 
general. Furthermore, the multiplicative inverse does not generally exists 
in this matrix algebra. However, one can easily obtain the inverse of any 
diagonal square matrix A with entries not equal to £ only on the diagonal. 
It is clear that with 



A 



I ai 



A 



-1 



e \ 



where > £ for all f = 1, . . . , ra, we have A (g) A~^ = A~^ ® A = E . 

Finally, the properties of matrix multiplication allow us to exploit the 
symbol A^ with a square matrix A and nonnegative integer p, as used in 
the conventional algebra: 



A^ = E, 



AP = AP 



-1 



^A = A^ 



A for p > 1 . 



p times 



2.2 A Linear Algebraic Equation 

Let us now examine a vector equation which will be encountered below in 
algebraic representations of tandem queueing system dynamics. For given 
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(n X n) -matrix A and column n-vector b, we consider the implicit equation 
in the n-vector x 

X = Ai^xQb, (1) 

which is generally identified analogous to the conventional linear algebra as 
a linear equation. 

The next lemma offers the solution of ([1]) for a particular class of matrices 
A . A detailed investigation of this and other linear equations in the general 
case can be found in [9|. 

Lemma 1. // there exists a nonnegative integer p such that AP = E , then 
([1]) has the unique solution 

p-i 

x = Y,e^'^b. (2) 

i=0 

Proof. Recurrent substitutions of x from ([T|) into the right-hand side of ([T]) 
and trivial algebraic manipulations give 

x = A^x®b = A0{A0x®b)®b = A'^0x®{E®A)^b 

= ■■■ = AP (S)xe{E(BAe---e a^^'^) <E) b. 

With the condition A^ = £ , we immediately arrive at (l2|) . It is also evident 
from the above calculations that the obtained solution is unique. □ 

As examples of the matrix A satisfying the condition of Lemma 1, we 
consider either lower or upper triangular matrices which have, in addition, 
the entries equal to e on the main diagonal. Specifically, it is not difficult 
to verify that for the matrix A with entries 

{Oj > e, ii i = j + 1 
e, otherwise, 

which will appear in the next sections, it holds that A^ ^ £ , for p = 
1, . . . ,n - 1, and ^" = 

3 Representation of Tandem System Dynamics 

We start with the scalar max-algebra equation representing the dynamics 
of a single-server queue, and then extend it to vector equations associated 
with tandem systems of queues. The queue is assumed to have a buffer with 
infinite capacity, and to operate under the first-come, first-served queue 
discipline. 

We denote the k^^ arrival epoch to the queue and the k^^ departure 
epoch from the queue by a{k) and d{k) respectively. The service time 
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of customer k is represented by ■ Under the conditions that the queue 
starts operating at time zero and it has no customers at the initial time, the 
dynamics of the system can be readily described in terms of a{k) and d{k) 
as state variables, by the following ordinary algebraic equation [UllllG]: 

d{k) = ma,x(a{k) , d{k — 1)) + r^. 

By replacing the usual operation symbols by those of max-algebra, one can 
rewrite this equation in its equivalent form as 

d{k) = Tk(E) a{k) (g) d{k - 1). 

Let us now suppose that there is a system of n single-server queues with 
infinite buffers. Furthermore, let ai{k) and di{k) be the k^^ arrival and 
departure epochs, and Tik be the service time of the k^^ customer for queue 
i. With the vector-matrix notations 



a{k) 



/ ai{k) 
\ an{k) 





/ d,{k) \ 




e \ 


d{k) = 










\ dn{k) I 


\ ^ 


Tnk J 


d{0) = (e, . . . , 


e)-^ , and d{k) = (e, . . 




for all /c < 



we may represent the dynamics of the whole system by the vector equation 

d{k) = Tk(E)a{k)®Tk(S>d{k-l). (3) 

Note that equation ^ is quite general; it may be appropriate for a 
variety of single-server queueing systems with infinite buffers, and not just 
for tandem systems. We will show how this equation can be refined in the 
case of two particular tandem systems in the next subsections. 

3.1 Closed Tandem Systems 

Consider a closed system of n single-server queues in tandem, with infinite 
buffers. In the system, the customers have to pass through the queues 
consecutively so as to receive service at each server. After their service 
completion at the n^^ server, the customers return to the first queue for a 
new cycle of service. It is assumed that the transition of customers from a 
queue to the next one requires no time. 

There are a finite number of customers circulating through the system. 
We suppose that at the initial time, all servers are free, whereas the buffer 
at the i^^ server contains Cj customers, i = 1, . . . ,n. It is easy to see that 
in this case, the arrival times of customers can be defined as 

a.(fc) = |^"^^-^^)' '^' = ' (4) 
\di-i{k-Ci), ifi = 2,...,n. 
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for all /c = 1,2,.... 

Let us now assume for simplicity that q = 1 for all i = 1, 
this assumption, one can rewrite (jH) in the vector form 



, n. With 



a{k) = F d{k - 1), where F 



e 



e 
e 



£ 



\ e e £ j 

Substitution of this expression for a(A;) into ([3]) leads to 

d{k) = Tk®F ®d{k-\)®Tk®d{k-\)=%®{F ®F)® d(k - 1). 

We may now describe the dynamics of the closed system under consideration 
by the linear state equation 



d{k) = Tfc (g) d(k - 1) 
with the state transition matrix 

Tk = Tk®{F® E) 



(5) 



£ 
£ 



Tlk \ 

£ 



£ 



Tn-lk £ 
"Tnk 'Tnk / 



(6) 



One can see that equation dSJ offers a very convenient way of calculating 
successive state vectors d{k) in the above system, by performing simple 
algebraic operations. Moreover, it is not difficult to understand that ^ can 
be readily extended to represent closed tandem systems with arbitrary Cj, 
< Q < cxD . A usual technique based on the employment of a modified state 
vector d{k) which integrates several consecutive state vectors of the original 
model, can be exploited to get ([5]) as a representation of such systems. 

To illustrate, let us suppose that at the initial time, there are q = 2 
customers in the buffer at each queue i, i = 1, . . . ,n. In this case, starting 
from (j3D, we may define the vector of the k^^ arrival epochs as 



and then rewrite 



a{k) = F(g)d{k-2) 
in the form 



dik) =Tk®d{k-l)®Tk(^F(^ d{k - 2) 



d{k) 



Finally, with the new state vector d{h) = 
to arrive at the state equation 

d{k) =fk(S>d{k-l 



d{k-l) 



, it is not difficult 
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with the state transition matrix 

^'=[ E £ 



3.2 Open Tandem Systems 



The purpose of this subsection is to extend equation ([5|) to the representation 
of open tandem queueing system models. In a series of single-server queues 
with infinite buffers, let us assign the first queue for representing an external 
arrival stream of customers. Each customer that arrives into the system has 
to pass through queues 2 to n, and then leaves the system. 

As in the closed system above, we suppose that Ttk represents the k^^ 
service time at queue i, i = 2, . . . ,n, whereas Tik now denotes the inter- 
arrival time between the k^^ customer and his predecessor in the external 
arrival stream. At the initial time, all servers are assumed to be free of 
customers, and, except for the first server, their buffers are empty; that is, 
Cj = 0, i = 2, . . . ,n. Finally, we put ci = oo to provide the model with the 
infinite arrival stream. 

It is clear that we have to define the arrival epochs in the system as 



ai{k) 



e, if i = 1 

di^i{k), if i = 2, 



(7) 



for all A; = 1, 2, ... . Proceeding to vector notation, we get 

/ e 

a{k) = G®d{k), where G = 

\ e e e J 

With this vector representation, equation Q takes the form 

d{k) = Tk®G^ d{k) ®Tk®d{k-l). 

The above equation can be considered as an implicit equation in d{k) , 
having the form of ([1]) with A = Tk®G and b = Tk® d{k — 1) . Moreover, 
one can readily see that the matrix 



Tk®G 



( e 

T2k 



\ 



I 



looks just like that presented in Subsection 2.2 as an example; therefore, it 
satisfies the condition of Lemma 1. 
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By applying Lemma 1, we obtain the solution 



n-1 



i=0 



Clearly, we have arrived at equation ([5]), where the transition matrix is 
defined as 



n-1 



i=0 

It is not difficult to verify that [11] 
Tk 

( T\k e e ... e \ 

T2k ^Tik T2k £ ... £ 

Tn-lfc (8) • • • (8) Tifc Tn-lk ® ' ' ' ^ T2k T^-lk ■ ■ ■ ® T^k £ 

r„fc ® • • • 18) Tifc r„fc (8 • • • 18 T2k Tnfc <8 • • • (8 Tsfc ... r„fc / 



(8) 



4 Tandem Queues with Finite Buffers and Block- 
ing 

Suppose now that the buffers of the servers in the open tandem system de- 
scribed above have limited capacity. Consequently, servers may be blocked 
according to one of the blocking rules. In this paper, we restrict our consid- 
eration to manufacturing blocking and communication blocking, which are 
most commonly encountered in practice [U [2] . 

Consider an open system with n servers in tandem, and assume the 
buffer at the i^^ server, i = 2, ...,n, to have capacity 6j, < 6j < cxd. 
We suppose that the buffer of the first server, which is the input buffer of 
the system, is infinite. Below is shown how the dynamics of tandem systems 
which operate according to both manufacturing and communication blocking 
rules, can be described by state equation ([5]). 

4.1 Manufacturing Blocking 

Let us first suppose that the dynamics of the system follows the manufac- 
turing blocking rule. Under this type of blocking, if upon completion of a 
service, the i^^ server sees the buffer of the {i + 1)*^* server full, the former 
server cannot be freed and has to remain busy until the (i + 1)*^* server 
completes its current service to provide a free space in its buffer. Clearly, 
since the customers leave the system upon their service completion at the 
j^th ggj-ygj. |;]2is server cannot be blocked. 
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It is not difficult to understand that the dynamics can be described by 
the ordinary scalar equations [H El [6] 

di{k) = max(max(ai(A;), di{k - 1)) + Tj^, di+i{k - bi+i - 1)), 

i = 1, . . . , n — 1, 
dn{k) = max{an{k),dnik - 1)) + Tnk, 

where ai{k) , i = 1, . . . ,n, are still defined by d?]). With max-algebra, one 
can readily rewrite these equations as 

diik) = ai{k) Tik di{k - 1) (B di+i{k - bi+i - 1), 
i = 1, . . . ,n — 1, 

dn{k) = Tnk <^ a„(/c) Tnk ^ dn{k - 1). 

Assuming 6j = 0, i = 2, . . . , n, for simplicity, we get the above set of 
equations in the vector form 

d{k) =Tk0 a{k) ®Tk®d{k-l)®G'^ <^ d{k - 1) 

= Tk(^ a{k) {Tk e G^) «) d{k - 1), 

where denotes the transpose of the above introduced matrix G . With 
a{h) = G (8) d{k) , we have the equation 

d{k) = Tk®G® d{k) e {Tk G^) ® d{k - 1). 

As in the case of the open tandem system with infinite buffers, we may 
apply Lemma 1 to solve this equation for d{k) , and obtain 

n-1 

d{k) = Y,e('Tk^Gy®{G^®Tk)(®d{k-l), 

i=0 

which equals ^ with the transition matrix 

n-1 

Tk = Y.^{Tk(^Gy®{G^(BTk). 

i=0 

Calculation of T^, [TT] leads us to 
Tk 

I Tik e e ... e \ 

T2fc (8) Tik T2k e e 

Tn-lk ®---(^Tik Tn-lk ® ' ' ' ® T2k T^-lk ® ' ' ' ® T^k e 
\ Tnk® ■■■® Tik Tnk®---®T2k Tnk ® ' ' ' ® T^k ■■■ Tnk ) 

(9) 
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Note that the matrices in ([8]) and ([9]) differ only in the elements of the upper 
diagonal adjacent to the main diagonal, which become equal to e in Q. 

Let us now derive equation 1^ for the system with the capacity of each 
buffer bi = 1, i = 2, . . . ,n. Clearly, the dynamics of the system can be 
described by the equation 

d{k) = Tk®G® d{k) ®Tk® d{k - 1) O d{k - 2). 

The application of Lemma 1 leads us to the equation 

n-l 

d{k) = ^eCTfc ^ Gy O (Tfc ® d{k -1)(BG^® d{k - 2)). 

i=0 

Finally, with the notations 

m = ( ) and Sk = J^eiTk G)\ 

we arrive at ([5]), taking the form 

dik) = l ^ ^ j(®d{k-l). 

4.2 Communication Blocking 

We now turn to a brief discussion of systems operating under the communi- 
cation blocking rule. This type of blocking requires a server not to initiate 
service of a customer if the buffer of the next server is full. In this case, 
the server remains unavailable until the current service at the next server is 
completed. 

Let us suppose that the above system follows communication blocking. 
The dynamics of this system can be described by the equations [U [B] 

di{k) = max{ai{k),di{k - l),di+i{k - bi+i - 1)) + nk, 

i = 1, . . . , n — 1, 
dn{k) = max(a„(A;), d„(A; - 1)) + T„fc, 

or, equivalently, by the max-algebra equations 

diik) = Tik (8) ai{k) Tik <^ di{k - 1) ® Tik di+i{k - bi+i - 1), 

i = 1, . . . , n — 1, 
dn{k) = Tnk an{k) T„fc (i„(A; - 1). 

For a particular system with 6j = 0, i = 2, . . . , n, we can write in the 
same manner as for manufacturing blocking 

d{k) =Tk<S) a{k) ® Tk ® d{k - 1) ® Tk <S) G'^ (S) d{k - 1) 

= Tk®G® d{k) ®Tk<^{E® G^) d{k - 1). 
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Furthermore, we have the solution 



n-l 



dik) = Y^^iTk GY ®Tk<^{E® G^) ® d{k - 1) 



i=0 



which takes the form of ([5]) with the matrix 



n-l 



Tk = J2^{Tk GY ^3 Tfc e G^) 




T2k 



e 



e 



e 



\ 



Tn-lk (8) • • • ® Tifc Tn-lk <^ ■ ■ ■ Tik T^-lfc • • • T2fc ... r„_ifc 

\ r„fc (g) • • • (g) rife r„fc (g) • • • (g) rife r„fc (g • • • (g r2fc ... r„fe ® r„_ifc / 



Let us finahy consider the system with the capacity of each buffer bi = 1, 
i = 2, . . . ,n. As for manufacturing blocking, we may represent the dynamics 
of the system by the equation 



5 Representation of Tandem System Performance 

In this section, we show how performance measures of open tandem systems 
may be represented based on the max-algebra models. We consider the 
measures representing the system time and the waiting time of customers, 
and give the corresponding linear max-algebra equations which allow us to 
describe and calculate these criteria in a simple way. The ordinary repre- 
sentation of other measures which one normally chooses in the analysis of 
queueing systems, including the utilization of a server, the number of cus- 
tomers at a queue, and the queue length, involves the operation of division 
[UElllllS], and they therefore cannot be expressed through linear equations 
in max-algebra. 

5.1 The System Time of Customers 

For the open tandem system with infinite buffers described above, let us 
define the vector of the system (sojourn) times of the k^^ customer as s{k) = 
{si{k), . . . , Sn{k)Y' , where Si{k) denotes the time required for the customer 
to pass through all queues up to and including queue i. Since the first 
queue represents the external arrival stream of customers, we do not have 



(10) 
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to include the time spent in this queue, and we therefore have si{k) = 0. 
Furthermore, the above definition of the system time leads to the equations 

Si{k) = di{k) - di{k), i = l,...,n. 

With max- algebra vector notations, we may rewrite these equations as 

s{k) =d{k)0d^\k). (11) 

It follows from ([5]) and the equality di{k) = di{k — 1) ^ Tik that 

s{k) =Tk® d{k -I)® d^^{k - 1) t~,}. 

By applying (jlip with k replaced by A; — 1 , we finally have 

s{k) = Uk(E)s{k-l), (12) 

where Uk = t'j} (g) Tk ■ 

As it is easy to see, the above dynamic equation is appropriate for rep- 
resenting the system time in open tandem systems with both infinite and 
finite buffers. For particular systems, equation (|12p will differ only in the 
matrix T^, and thus in the matrix Uk, inherent in their associated dynamic 
state equations. Finally, note that a lower triangular transition matrix Tk 
will result in a matrix Uk of the same kind. One can consider the open 
tandem system with infinite buffers as an example. 



5.2 The Waiting Time of Customers 

Consider the open tandem system with infinite buffers once again. The 
system time of a customer in this system consists of his service time and the 
waiting time. Therefore, introducing the symbol Wi{k) for the k^^ customer 
to denote the total time spent on waiting for service at the queues from 1 
to i , we have the equations 

si{k) = wi{k) = 0, 

i 

Si{k) = Wi{k) + '^Tjk, i = 2,3, ...,n. 

In order to represent the relation between the system and waiting times 
in max-algebra vector form, let us further introduce the vector w^k) = 
{wi{k), . . . ,Wn{k)y , and the diagonal matrix 



Pk 



( rik £ 

e Tlfc ® T2k 



\ 
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Clearly, we may now write 

s{k) = (g) Pk®w{k). 

Since the multiplicative inverse exists for the matrix , it is not difficult 
to resolve the above equation for w{k) : 

w{k) = Tik®Pk^ (S> s{k). 

With (fT2]l . we successively obtain 

w{k) = rife (g) P^^ ® (g) Tfe O s{k - 1) 

= P-^ (g) Tfe O T-^_^ (g) Pfe_i (g) (ru._i ® p-\ O s{k - 1)) 

= Tffcii O (g Tfe (g) Pfe_i (g) - 1). 

Finally, we arrive at the dynamic equation 

w{k) = Vk,k-i0w{k-l), (13) 

with the transition matrix Vfe,fe_i = ® P]^^ (g) ?fe (g Pfe_i . 

Note that, in general, equation (fT^ can be extended to tandem queues 
with finite bufi^ers and blocking. In that case, however, the quantities Wi{k) 
will include not only the time spent on waiting for service, but also the 
blocking time of servers. 

6 Serial and Parallel Simulation of Tandem Queues 

In this section, we briefly discuss serial and parallel simulation algorithms 
based on the algebraic models presented above, and outline their perfor- 
mance. We take, as the starting point, the state equation 

d{k) =Tk(S)d{k-l). (14) 

with the matrix Tfe defined by ([8]). Clearly, calculations with other matrices, 
including performance evaluation via equations ()12p and ()13p . will normally 
differ little in complexity. 

Furthermore, we assume as in [H El [7] that the service times Tik , which 
are normally defined in queueing system simulation as realizations of given 
random variables, are available for all i = 1, . . . , n and k = 1,2 ... , when 
required. We will therefore concentrate only on procedures of evaluating the 
system state vectors, based on ([HI)- 
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6.1 Simulation with a Scalar Processor 



Let us first consider time and memory requirements when only one scalar 
processor is available for simulation of tandem systems. The related serial 
algorithm consists of consecutive steps, evaluating new state vectors. The 
k^^ step involves determination of the transition matrix and multipli- 
cation of this matrix by the vector d(k — 1) to produce the vector d(k) . 
However, particular algorithms designed for computation with matrices 
dH), (HJ, and (fTOj) may execute the step in different ways, according to the 
structure of the matrices, so as to reduce time and memory costs. 

Assume that a system is simulated until the K^^ service completion at 
server n, and denote the overall number of the operations © and ® to be 
performed within the k^^ step to evaluate the matrix T^, and to compute the 
product Tfc (g) d{k — 1) respectively by A''i and In that case, the entire 
simulation algorithm will require = K{Ni + N2) operations, ignoring 
index manipulations. Finally, we denote by M the number of memory 
locations involved in the computations. 

Let us now consider the open tandem system with infinite buffers, and 

(k) 

denote the entries of its transition matrix T^. by t^j , i = j = 

1, . . . ,n. Taking into account the lower triangular form of the matrix 
defined by ([8]), a serial algorithm for calculating K successive vectors d(k) 
may be readily designed as follows. 

Algorithm 1. 

For i = 1, . . . ,n, do dj(0) < — e. 
For k = l,...,K, do 
for i = 1, . . . , n, do 

t^'^ ^ t;- 

for J = 1, . . . , i — 1, do 

for i = 1, . . . , n, do 

di{k) ^ t^^ ® di{k - I)- 
for j = 2, . . . ,i, do 

di{k) ^ di{k) © 0dj{k-l). 

It is not difficult to see that for the algorithm, A^i = n{n+l)/2 and N2 = 
r? operations are required; it entails M = n(n + 5)/2 memory locations. 
With the total number of operations A^ = 0(n(3n + l)Ar/2) involved in 
calculating K state vectors, the algorithm proves to have a performance 
comparable to other serial simulation procedures [UIZ!- As an example, the 
serial algorithm in [7j, designed for simulation of open tandem systems with 
infinite buffers can be considered. This algorithm allows one to compute the 
AT**^ departure time at queue n in 0{2K{n + 1)) operations; however, it is 
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actually a scalar algorithm intended to obtain only the value of dn(k) rather 
than the whole vector d{k) . Therefore, one has to multiply its performance 
criterion by n so as to provide a proper basis for a comparison between 
algorithms. It should be noted in conclusion that calculations with matrix 
([6]) can be performed using only N = 0{2Kn) operations and M = 0(3n) 
memory locations. 

6.2 Simulation with a Vector Processor 

Suppose now that simulation is executed on a vector processor equipped with 
vector registers of a length large enough for the processing of n -vectors. 
With the notations = (t -^^ , • • • , ) , i = 1 , . . . , n , we may write the 
following vector modification of Algorithm 1: 

Algorithm 2. 



First note that implementation of the vector processor makes it possible 
to compute matrix ([8]) by performing only Ni = n vector operations. Fur- 
thermore, evaluation of each element of the vector d{k) from ()14p actually 
involves both componentwise addition of a row of the matrix to the vec- 
tor d{k — 1) , and determination of the maximum over the elements of this 
vector sum. 

With the vector processor, one can add two vectors together in a single 
operation. It follows from the triangular form of ([8]) that to compute the i^^ 
element of d{k) , one actually has to perform no more than i maximizations. 
By applying the recursive doubling method [ISJ , we may obtain the maximum 
over i consecutive elements of a vector in log2 i operations. For all entries of 
the vector d{k) , we therefore get A''2 = n+log2 1 + - • •+log2 n = n+log2(n!) . 

Finally, evaluation of K state vectors requires N = 0{K{\og2{'n\) + 2n)) 
operations of a vector processor. As is easy to see, vector computations allow 
us to achieve the speedup = 0(n(3n + l)/(log2(n!)/2 + n)) in relation to 
the sequential procedure discussed in the previous subsection. 

6.3 Parallel Simulation of Tandem Systems 

We conclude this section with the discussion of a parallel simulation al- 
gorithm intended for implementation on single instruction, multiple data 



d(0) ^ (e,...,e)'^. 
For = 1, . . . ,K, do 

Tk < — £; 

for i = 1, . . . ,n, do 
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(SIMD) parallel processors. Other parallel algorithms based on the alge- 
braic models can be found in [2j. 

Let P be the number of processors available, and P >n. The algorithm 
consists of L = [X/P] steps, where \x\ denotes the smallest integer greater 
than or equal to x, and it can be described as follows. 

Algorithm 3. 

d(0) ^ (e,...,er. 
for / = 1, . . . , L, do 

in parallel, for i = (/ - 1)P + 1, (/ - 1)P + 2, . . . , min(/P, K), do 
evaluate Xi; 

for i = (/-l)P + l,(/-l)P + 2, ...,min(/P,E:), do 

in parallel, for j = 1, . . . , n, do dj{i) < — t^*"* (8) d{i — 1). 

To estimate the performance of this algorithm, first note that each step 
starts with parallel evaluation of P consecutive transition matrices, which 
entails A^i = n(n + l)/2 parallel operations. Next follows the determination 
of P consecutive state vectors, with each vector evaluated in parallel by 
computing its elements on separate processors. Since evaluation of a vector 
element involves n ordinary additions and the same number of maximiza- 
tions, one has to perform 2n parallel operations so as to determine the whole 
vector. Therefore, we have N2 = 2Pn operations required for computing P 
vectors. 

Finally, the entire algorithm entails N = L{n{n + l)/2 + 2Pn) par- 
allel operations. It is not difficult to understand that with P > n, we 
get N = 0{3Kn/2) as ^ cxo. Moreover, by comparison with the 
performance of the above sequential procedure, one can conclude that for 
P = n and K sufficiently large, the parallel algorithm achieves the speedup 
Sp = 0(3P/5). 

7 Conclusions 

The max-algebra approach provides a very convenient way of representing 
tandem queueing system models. The algebraic models describe the dy- 
namics of the systems through linear max-algebra equations similar to those 
being studied in the conventional linear systems theory. The approach there- 
fore offers the potential for extending classical results of both conventional 
linear algebra and systems theory to analyze queueing systems. Moreover, it 
provides the basis for applying methods and algorithms of numerical algebra 
to the development of efficient procedures for the queueing systems simula- 
tion, including the algorithms designed for implementation on parallel and 
vector processors. 
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Since the closed form representation of the dynamics of more comphcated 
queueing systems including the G/G/m queue and queueing networks, nor- 
mally involves three operations, namely addition, maximization, and min- 
imization (see, e.g., [3 [6]), these systems cannot be described in terms of 
max-algebra. In that case, however, appropriate algebraic representations 
can be obtained using minimax algebra [HI llOj. 
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